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A  NONRESONANT  TECHNIQUE  TO  ESTIMATETHE  MECHANICAL 
MATERIAL  PROPERTIES  OF  A  VISCOELASTIC  CYLINDER 


1.  INTRODUCTION 

Measuring  the  Young’s  modulus,  shear  modulus,  and  Poisson’s  ratio  of  materials  is 
important  because  these  parameters  contribute  significantly  to  the  static  and  dynamic  response  of 
a  structure. 

Resonant  techniques  [1-7],  which  have  been  used  to  identify  and  measure  moduli  for  many 
years,  are  based  on  comparing  the  measured  eigenvalues  of  a  structure  to  predicted  eigenvalues 
from  a  model  of  the  same  structure.  One  limitation  is  that  the  model  of  the  structure  must  have 
well-defined  (typically  closed-form)  eigenvalues.  A  second  concern  is  that  resonant  techniques 
allow  measurements  only  at  natural  frequencies. 

Comparison  of  analytical  models  to  measured  frequency  response  functions  [8-12]  is  another 
approach  for  estimating  the  mechanical  properties  (stiffness  and  loss)  of  a  structure.  When  the 
analytical  model  agrees  with  one  or  more  frequency  response  functions,  the  parameters  used  to 
calculate  the  analytical  model  are  considered  accurate.  However,  if  the  analytical  model  is 
formulated  with  a  numerical  method,  a  comparison  of  the  model  to  the  data  can  be  difficult  due 
to  the  dispersion  properties  of  the  materials. 

Other  methods  that  measure  Young’s  modulus  require  strain  gauges  [13-14]  to  be  affixed  to 
the  cylindrical  test  specimen.  However,  the  mounting  of  such  gauges,  which  normally  requires 
that  the  gauge  be  glued  to  the  specimen,  stiffens  the  material  locally,  causing  an  especially 
adverse  impact  on  the  estimate  of  stiffness  and  loss  for  soft  viscoelastic  materials. 

Still  another  technique  for  measuring  stiffness  and  loss  is  to  deform  the  material  and  measure 
the  resistance  to  the  indentation  [15].  The  drawback  here  is  that  the  specimen  can  be  physically 
damaged  if  the  deformation  causes  the  sample  to  enter  the  plastic  region  of  deformation. 
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This  report  describes  the  development  of  a  nonresonant  inverse  method  for  measuring  the 
complex  Young’s  modulus  and  complex  shear  modulus  of  a  cylindrical  test  specimen  containing 
a  mass  on  its  end.  The  technique  is  best  used  when  the  end  mass  will  change  the  dynamic 
response  of  the  system,  which  typically  occurs  when  an  end  accelerometer  that  measures 
response  is  attached  to  the  test  specimen.  (It  should  be  noted  that  the  larger  the  mass,  the  more 
significant  will  be  the  response  of  the  specimen.) 

In  the  next  sections,  the  linear  and  angular  equations  of  motion  for  the  test  specimen  with  the 
mass  are  derived.  The  inverse  method  allows  two  sets  of  experimental  measurements  in  the  form 
of  transfer  functions  to  be  recorded  for  both  the  linear  and  rotational  configurations.  These  sets 
are  combined  to  yield  closed-form  values  of  Young’s  modulus  and  shear  modulus  at  each 
frequency  for  which  a  measurement  is  made.  The  calculation  of  Poisson’s  ratio,  based  on  a 
material  constitutive  equation,  can  then  be  made  from  these  previous  measurements  of  moduli. 
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2.  MEASUREMENT  OF  YOUNG’S  MODULUS 


Measuring  the  Young’s  modulus  of  the  test  specimen  begins  by  mounting  one  end  of  the 
specimen  to  a  mechanical  shaker  and  attaching  the  other  end  to  a  mass.  Next,  one  accelerometer 
is  attached  to  the  shaker  head  and  a  second  is  attached  to  the  end-mounted  mass,  with  the 
measurement  axis  of  both  accelerometers  in  the  x  direction.  The  laboratory  configuration  for  this 
experiment  is  shown  in  figure  1.  As  can  be  seen  from  the  figure,  the  shaker  inputs  energy  into  the 
test  specimen  in  the  form  of  a  linear  translation  that  initiates  a  compressional  wave.  The  speed 
and  loss  of  this  wave  is  measured  with  the  two  accelerometers,  and  Young’s  modulus  is  then 
calculated  from  the  result.  This  measurement  (estimation)  process  is  further  described  below. 

The  system  model  represents  the  test  specimen  (bar)  attached  to  the  mechanical  shaker  at 
x  =  0  and  the  mass  at  x  =  L,  with  the  two  attached  accelerometers  measuring  the  acceleration 
levels  at  each  end  of  the  bar.  The  linear,  second-order  wave  equation  modeling  displacement  in 
the  bar  is 

d2u(x,t)  Ed2u(x,t)_ 

dt2  P  dx2  '  .  W 


where  u(x,t)  is  the  axial  particle  displacement  at  location  x  (m)  and  time  t  (s),  p  is  the  density  of 
the  bar  (kg/mJ),  and  E  is  the  unknown  frequency-dependent,  complex  Young’s  modulus  of 
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elasticity  (N/m~ ),  which  is  yet  to  be  determined.  The  boundary  at  x  =  0  is  modeled  as  a  fixed 
end  with  harmonic  motion  and  is  expressed  as 

m(0,0  =  f/0exp(ift»)  ,  (2) 


where  ©is  the  frequency  of  excitation  (rad/s),  Uq  is  the  amplitude  (m),  and  i  is  the  square  root 

of  -1.  The  boundary  at  x  =  L  is  formulated  by  matching  the  force  at  the  end  of  the  bar  to  the  end 
mass  multiplied  by  its  linear  acceleration  and  is  expressed  as 


AE 


du(L,t) 

dx 


-m 


d2u(L,t ) 
dt2 


(3) 
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Accelerometer 


Figure  1.  Laboratory  Configuration  for  Measuring  Young’s  Modulus 


where  A  is  the  cross-sectional  area  of  the  bar  (nr )  and  m  is  the  mass  at  the  end  of  the  bar  (kg). 


Equation  (1)  can  be  rewritten  in  the  spatial  domain  as 


dzU(x,o) )  ,2tt,  x  ^ 

- - —  +  kEU(x,(o)  =  0  , 

cbr 


(4) 


where  U(x,o) )  is  the  temporal  Fourier  transform  of  the  axial  displacement  and  kE  is  the  complex 
compressional  wavenumber  (rad/m),  which  is  equal  to 


kE  =(q/t]e/p  . 


(5) 


Similarly,  equation  (2)  becomes 


U(0,Q))  =  Uq  , 


(6) 


and  equation  (3)  becomes 


.„d  U(L,C0)  2  ri,r  \ 

AE - =  mco  U(L,a)) 

dx 


(7) 


The  solution  to  equation  (4)  is 


U(x,Q) )  =  R(co)cos(kEx)  +  S(co)sm(kEx)  , 


(8) 


where  R  and  S  are  wave  propagation  constants.  Applying  boundary  conditions  (6)  and  (7)  to 
equation  (8)  and  then  writing  the  solution  as  a  transfer  function  in  the  form  of  a  ratio  between  the 
displacement  at  both  ends  produces 


U(L,o) ) 


U, 


0 


1 


cos(kEL)~  p(kEL)sin(kEL) 


(9) 
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where  p  is  the  ratio  of  the  attached  mass  to  the  bar  mass  and  is  equal  to 


m 


where  M  is  the  mass  of  the  bar  (kg)  expressed  as 


(10) 


M  =  pAL  . 


(ID 


The  transfer  function  in  equation  (9)  represents  experimental  data  and  is  a  function  of 
unknown  wavenumber  kE .  The  inversion  of  the  two  transfer  functions  obtained  from  the 
different  attached  masses  will  allow  both  sets  of  data  to  be  combined,  yielding  a  closed-form 
solution  of  kE  and  then  £  as  a  function  of  (0.  The  theoretical  forms  of  these  transfer  functions  are 


cos  (kEL)  -  Pi(kEL)sin(kEL) 

and 


U(L,cq) 


Ut 


0 


=  2i(o»  = 


cos  (kEL)  -  p2  (kEL)  sin  (kEL)  J 

where  the  subscript  1  denotes  the  first  attached  mass  and  the  subscript  2  denotes  the  second 
attached  mass.  Although  the  theoretical  forms  given  in  equations  (12)  and  (13)  are  end 
displacement  divided  by  input  displacement,  the  transfer  function  created  by  end  acceleration 
divided  by  input  acceleration  (obtained  from  the  accelerometers)  is  identical  in  the  frequency 
domain.  Writing  equations  (12)  and  (13)  as  a  function  of  (kEL)sin(kEL)  and  then  equating  them 
yields 


U(L,co) 


U{ 


0 


=  T2(co)  = 


cos  (kEL)  = 


t2^2  ~  T\fh 
W2^2~lx l) 


=  0  > 


(14) 
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where  0  is  a  complex  quantity.  The  inversion  of  equation  (14)  allows  the  complex  wavenumber 
to  be  solved  as  a  function  of  (j).  The  solution  to  the  real  part  of  kE  is 


— Arccos(s)  +  — 
2  L  2  L 


Re(^)  = 


—Arc  cos  (-^)  +  — 
2L  2  L 


n  even 


n  odd 


(15) 


where 


5  =  [Re(0)]2  +  [lm(0)]2 
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[Re(0)]2  +  [lm(^)]2 


}2-{2[Re(0)]2-2[Im(^))2-l}, 


(16) 


n  is  a  non-negative  integer,  and  capital  A  denotes  the  principal  value  of  the  inverse  cosine 
function.  The  value  of  n  is  determined  from  the  function  s,  which  is  a  cosine  function  with 
respect  to  frequency.  At  zero  frequency,  n  is  0.  Every  time  s  cycles  through  tt  radians  (180°), 
n  is  increased  by  1.  After  the  solution  to  the  real  part  of  kE  is  found,  the  solution  to  its 
imaginary  part  is  written  as 

Im  (kE)  =  -loge{ - - - 

L  e[cos[R  e(kE)L]  sin[R  e(kE)L] 

Once  the  real  and  imaginary  parts  of  wavenumber  kE  are  known,  the  complex- valued  modulus 
of  elasticity  can  be  determined  at  each  frequency  with 

2 

E((0 )  =  Re[£(to)]  +  iIm[£((w)]  =  - ^ - —  .  (18) 

'  [Re(£  £ )  +  i  Im(  A:  E)]2 

Equations  (12)-(18)  produce  an  estimate  of  Young’s  modulus  at  every  frequency  for  which  a 
measurement  is  made. 
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3.  NUMERICAL  SIMULATION  OF  YOUNG’S  MODULUS  MEASUREMENT 


Numerical  simulations  conducted  to  determine  the  effectiveness  of  this  method  use  the 
following  parameters  to  define  a  baseline  problem:  M  =  4.0  kg,  m\  =  0.4  kg,  m2  =  1.2  kg, 

L  =  0.254  m,  p  =  1200  kg/m3,  Re(£)  =  108  + 105/  N/m2 ,  and Im(£)  =  107  + 104/  N/m2 , 
where /is  frequency  in  Hertz.  From  these  values,  the  mass  ratios  /tj  and  /t2  are  computed  to  be 
0.1  and  0.3,  respectively. 

Figure  2,  which  corresponds  to  equations  (12)  and  (13),  is  a  plot  of  magnitude  (top)  and 
phase  angle  (bottom)  for  the  transfer  functions  7}  (ft))  and  r2(ft))  versus  frequency.  The  first 
transfer  function,  computed  with  an  attached  mass  of  0.4  kg,  is  depicted  with  x’s,  and  the  second, 
computed  with  an  attached  mass  of  1.2  kg,  is  shown  with  o’s. 

Figure  3,  corresponding  to  equation  (16),  is  a  plot  of  the  function  s  versus  frequency.  Note 
that  although  this  function  is  a  cosine  with  respect  to  frequency,  the  period  is  increasing  as 
frequency  increases.  Table  1  lists  the  values  of  n  versus  frequency,  which  were  determined  by 
inspection  of  the  figure.  Once  the  values  of  n  are  known,  the  modulus  values  of  E  can  be 
determined  from  equations  (15)-(18). 

Figure  4  shows  a  plot  of  real  and  imaginary  Young’s  modulus  values  versus  frequency.  The 
real  (actual)  values  used  to  make  the  transfer  functions  are  displayed  as  a  solid  line,  and  the  real 
(estimated)  values  are  displayed  as  x’s.  The  imaginary  (actual)  values  used  to  make  the  transfer 
functions  are  also  displayed  as  a  solid  line,  with  the  imaginary  (estimated)  values  displayed  as 
o’s.  The  estimated  values  agree  at  all  frequencies  with  the  actual  values.  This  result  is  expected 
because  there  is  no  noise  in  the  data  and  the  same  parameters  used  to  make  the  transfer  functions 
are  used  to  calculate  the  modulus  values  (i.e.,  no  error  is  introduced  when  calculating  the 
modulus  from  the  transfer  functions). 

Errors,  however,  are  introduced  into  this  method  during  the  measurement  of  the  transfer 
functions  and  system  parameters.  This  effect  can  be  simulated  by  adding  noise  to  the  data  with 
random  numbers  added  to  the  transfer  function,  as  shown  in 

Te (ft))  =  T(ft>)  +  e{ Re[T (ft))]<7 a  +ilm [T(co)ab]}  ,  (19) 
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Phase  Angle  (degrees) 


Magnitude 


where  e  is  the  amount  of  error  added  to  the  transfer  function  and  where  aa  and  <7£  are  random 
numbers  with  zero  mean  and  a  variance  of  one. 


Figure  5,  corresponding  to  equations  (12),  (13),  and  (19),  is  a  plot  of  the  transfer  functions 
T\((o)  and  ^(ft))  with  noise  versus  frequency.  The  value  of  e  used  to  create  this  figure  is  0.03, 
which  represents  approximately  3%  noise  in  the  data  set.  (If  this  test  is  properly  set  up  and 
executed,  a  3%  noise  level  would  be  considered  even  greater  than  the  actual  noise  levels  in  the 
experimental  data.) 


Figure  6  is  a  plot  of  the  function  s  with  and  without  additive  noise  versus  frequency.  The 
solid  line  is  the  function  with  e  set  equal  to  0.03,  and  the  o’s  connected  with  the  dotted  line  show 
the  function  without  noise  (also  seen  in  figure  3). 


Figure  7  plots  the  estimated  and  actual  values  of  Young’s  modulus  (E)  with  noise  versus 
frequency.  The  estimated  values  of  the  real  part  of  E  are  depicted  with  x’s,  and  the  actual  values 
of  the  real  part  of  E  are  shown  as  a  solid  line.  The  estimated  values  of  the  imaginary  part  of  E 
are  depicted  with  o’s,  and  the  actual  values  of  the  imaginary  part  of  E  are  shown  as  a  solid  line. 


The  effect  of  measurement  error  on  the  accuracy  of  the  calculation  of  E  can  also  be  studied 
with  Monte  Carlo  simulations.  Eleven  different  values  of  e  were  used  to  build  the  transfer 
functions  and  then  calculate  the  modulus  value.  Estimation  error  at  each  frequency  was  defined 
according  to 


Eact  )|  \^est  (fi^m  )|| 
max[£acf  (com 

)^est  (fom)]  ’ 


(20) 


where  a(com)  is  Young’s  modulus  estimation  error  at  the  mth  frequency,  the  subscript  act 
corresponds  to  the  actual  value,  and  the  subscript  est  corresponds  to  the  estimated  (calculated) 
value.  Once  this  value  is  known,  it  can  be  summed  across  M  frequencies  by 

,  M 

Pj  =  ,  (21) 

m= 1 
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where  fij  is  the  average  error  for  the  yth  Monte  Carlo  simulation.  Finally,  all  the  simulations  can 
be  summed  from 

1  J 

j  ,  (22) 

where  e  is  the  estimation  error  using  J  simulations.  Table  2  lists  transfer  function  error  (e)  versus 
Young’s  modulus  estimation  error  (e).  The  frequency  range  of  the  values  used  to  prepare  this 
table  is  2.5  to  5000  Hz,  the  number  of  frequencies  (M)  used  to  calculate  the  estimation  error  is 
100,  and  the  number  of  simulations  (7)  for  each  value  of  e  is  100. 

Figure  8  is  a  plot  of  Young’s  modulus  estimation  error  (e)  versus  frequency  for  a  single 
simulation  (7=1)  using  a  known  transfer  function  error  of  e  =  0.03.  Note  that  the  estimation 
error  for  Young’s  modulus  is  slightly  smaller  than  the  transfer  function  error  in  table  2.  Also 
note  that  the  larger  estimation  errors  tend  to  be  near  the  frequency  values  where  n  is  increasing. 
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Phase  Angle  (degrees) 


Frequency  (Hz) 


Figure  5.  Magnitude  and  Phase  Angle  for  Transfer  Functions  T\  and  T2 

With  Noise  Versus  Frequency 


Figure  7.  Real  and  Imaginary  Young’s  Modulus  Values 
Calculated  With  Noise  Versus  Frequency 


Table  2.  Transfer  Function  Error  Versus 
Young’s  Modulus  Estimation  Error 


Transfer 

Estimation 

Function  Error 

Error 

(e) 

(e) 

0 

0 

0.005 

0.004 

0.010 

0.008 

0.015 

0.013 

0.020 

0.019 

0.025 

0.023 

0.030 

0.026 

0.035 

0.031 

0.040 

0.034 

0.045 

0.038 

0.050 

0.041 

E  Estimation  Error 


4.  MEASUREMENT  OF  SHEAR  MODULUS 


Measuring  the  shear  modulus  of  the  cylindrical  test  specimen*  requires  mounting  the  shaker 
longitudinally  to  a  stiff  connecting  rod  (not  to  the  specimen),  which  is  connected  to  the  edge  of  a 
large  disk.  A  bearing  mounted  at  the  center  of  the  disk  allows  it  to  rotate  freely  about  its  center 
point.  The  test  specimen  is  rigidly  attached  to  the  middle  of  this  disk  in  such  a  manner  that  when 
the  disk  is  pushed  by  the  connecting  rod,  it  initiates  torsional  (or  rotational)  response  in  the 
specimen.  The  other  end  of  the  test  specimen  is  attached  to  a  second  disk  that  provides  rotary 
inertia  when  the  test  is  run.  Two  accelerometers  are  used  -  one  attached  to  the  first  disk  and  one 
attached  to  the  second  disk  -  with  the  measurement  axis  of  both  in  the  same  angular  direction  as 
the  disks.  Although  each  accelerometer  measures  translation,  this  value  can  be  converted  into 
angular  rotation  by  multiplying  the  recorded  value  by  the  distance  from  the  center  of  the  disk  to 
the  accelerometer.  The  laboratory  configuration  for  this  experiment  is  shown  in  figure  9. 

As  can  be  seen  in  the  figure,  the  shaker  inputs  energy  into  the  connecting  rod  in  the  form  of 
linear  translation.  When  the  rod  inputs  the  energy  into  the  disk,  the  disk  rotates  and  initiates  a 
shear  wave  in  the  test  specimen.  The  speed  and  loss  of  this  wave  can  be  measured  with  the  two 
accelerometers,  and  the  shear  modulus  can  be  calculated  from  the  result.  This  measurement 
(estimation)  process,  which  is  extremely  similar  to  the  measurement  of  Young’s  modulus,  is 
described  below. 


The  system  model  represents  a  cylindrical  test  specimen  (bar)  attached  via  a  disk  and 
connecting  rod  to  a  torsional  shaker  at  x  =  0  and  to  a  disk  with  rotary  inertia  at  x  =  L.  Each  disk 
has  an  accelerometer  to  measure  the  angular  acceleration  levels  at  its  end.  The  linear  second- 
order  wave  equation  modeling  angular  rotation  in  the  bar  is 


<?20(x,r)  G<?20(;r,O_o 

dt 2  p  dx2  ~ u 


(23) 


*As  with  Young’s  modulus,  two  sets  of  experimental  measurements  are  recorded,  with  the  disk  (at  x  =  L)  in  the 
second  set  having  a  different  rotary  inertial  value  than  the  one  in  the  first  set. 
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where  0(x,t)  is  the  angular  rotation  at  location  *  (m)  and  time  t  (s),  p  is  the  density  of  the  bar 
( kg/m  ),  and  G  is  the  unknown,  frequency-dependent,  complex  shear  modulus  of  elasticity 
(N/m  ),  which  is  yet  to  be  determined.  The  boundary  at  x  =  0  is  modeled  as  a  fixed  end  with 
harmonic  angular  motion  and  is  expressed  as 

0(0,0  =  0oexp(imO  ,  (24) 

where  (O  is  the  frequency  of  excitation  (rad/s),  0q  is  the  amplitude  (rad),  and  i  is  the  square  root 

of  -1.  The  boundary  at  x  =  L  is  formulated  by  matching  the  angular  force  (torque)  at  the  end  of 
the  bar  to  the  rotary  inertia  of  the  end  mass  multiplied  by  its  angular  acceleration  and  is 
expressed  as 

d0(L,t)  d2G(L,t) 

GIr^r  =  -J^-’  (25) 


where  Ip  is  the  polar  moment  of  inertia  in  the  cross-section  of  the  bar  (m4 )  and  J  is  the  rotary 
inertia  of  the  disk  at  the  end  of  the  bar  (kg«m2).  For  the  cylindrical  test  specimen,  the  polar 
moment  of  inertia  is 


(26) 


where  a  is  the  radius  of  the  specimen  (m).  For  a  cylindrical  disk,  the  rotary  inertia  is 


(27) 


where  r  is  the  radius  of  the  disk  (m)  and  m  is  the  mass  of  the  disk  (kg). 
Equation  (23)  can  be  rewritten  in  the  spatial  domain  as 
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(28) 


d  ®(x,a>)  ,  ,2 


dx 


+  kG@(x,co)  =  0  , 


where  0(x,  (0)  is  the  temporal  Fourier  transform  of  the  axial  displacement  and  kG  is  the  complex 
shear  wavenumber  (rad/m),  which  is  equal  to 


kG  =  (0/y[G/p  . 


(29) 


Similarly,  equation  (24)  becomes 


0(O,to)  =  0o  , 


(30) 


and  equation  (25)  becomes 


dQiUo)  =  /(y20(I  ?(y) 
F  dx 


(31) 


The  solution  to  equation  (28)  is 


0(x,tu)  =  X((0)cos(kGx )  +  Z(co)sm(kGx)  , 


(32) 


where  X  and  Z  are  wave  propagation  constants.  Applying  boundary  conditions  (30)  and  (3 1)  to 
equation  (32)  —  and  writing  the  solution  as  a  transfer  function  in  the  form  of  a  ratio  between  the 
rotation  at  both  ends  —  produces 


0(L,m) 


©o  [cos(kGD-MkGL)sin(kGL)j 


(33) 


where  A  is  equal  to 
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2/ 


(34) 


a  M 

and  M  is  the  mass  of  the  bar  (kg)  expressed  as 


M  =  pAL  .  (35) 

The  transfer  function  in  equation  (33)  represents  experimental  data  and  is  a  function  of 
unknown  wavenumber  kG .  The  inversion  of  the  two  transfer  functions  obtained  from  the 
different  attached  rotary  inertial  masses  will  allow  both  sets  of  data  to  be  combined,  yielding  a 
closed-form  solution  of  kG  and  then  G  as  a  function  of  to.  The  theoretical  forms  of  these  transfer 
functions  are 


0(L,to) 

©0 


=  $!(©)  = 


_ 1 _ 

cos  (kG  L)  -  A1  (kG  L)  sin{kG  L ) 


(36) 


and 


eo»_  c  i 

- —  So  ((D)  —  - 

©o  z  lcos(kG  L)  -  X2  ikG  L)  &in(kG  L) 


(37) 


where  the  subscript  1  denotes  the  first  attached  rotary  inertial  mass  and  the  subscript  2  denotes 
the  second  attached  rotary  inertial  mass.  Although  the  theoretical  forms  given  in  equations  (36) 
and  (37)  are  end  rotation  divided  by  input  rotation,  the  transfer  function  created  by  end  rotational 
acceleration  divided  by  input  rotational  acceleration  (which  is  obtained  from  the  accelerometers) 
is  identical  in  the  frequency  domain.  Writing  equations  (36)  and  (37)  as  a  function  of 
(kGL)sin(kGL)  and  then  equating  them  yields 


cos (kGL)  - 


^2^2  ~^1^1 
^1^2  (^2  ~^l) 


=  <P  , 


(38) 
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where  <p  is  a  complex  quantity.  The  inversion  of  equation  (38)  allows  the  complex  wavenumber 
to  be  solved  as  a  function  of  (p .  This  solution  to  the  real  part  of  kG  is 


Re(*G)  = 


1  .  .  .  mn 

— Arccos(r)  + -  m  even 

2  L  2  L 

1  .  .  .  mn  ,, 

— Arccos(-r)  + m  odd 

2  L  2  L 


(39) 


where 


r  =  [Re(<p)]2  +  [1mm2  ~  ^/{[Re(<p)]2  +[Im(<p)]2}2-{2[Re((p)]2  -2[lm((p)]2  -l}  , 

(40) 

m  is  a  non-negative  integer,  and  capital  A  denotes  the  principal  value  of  the  inverse  cosine 
function.  The  value  of  m  is  determined  from  the  function  r,  which  is  a  cosine  function  with 
respect  to  frequency.  At  zero  frequency,  m  is  0.  Every  time  s  cycles  through  7T  radians  (180°), 
m  is  increased  by  1.  After  the  solution  to  the  real  part  of  kG  is  found,  the  solution  to  its 
imaginary  part  is  written  as 


1  w  f  RefoO 


Im(<p)  ) 


Im(£G)  Llo§ejcos[Re(^)L]  sin[Re(£G)L]J  • 


(41) 


Once  the  real  and  imaginary  parts  of  wavenumber  kG  are  known,  the  complex-valued  modulus 
of  elasticity  can  be  determined  at  each  frequency  with 

2 

G(co)  =  Re[G(<o)]  +  i  Im[G(<o)]  = - ^ - -  .  (42) 

[Re(&G )  +  i  Im(fcG  )]2 

Equations  (23)-(42)  produce  an  estimate  of  shear  modulus  at  every  frequency  for  which  a 
measurement  is  made. 
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5.  NUMERICAL  SIMULATION  OF  SHEAR  MODULUS  MEASUREMENT 


Numerical  simulations  conducted  to  determine  the  effectiveness  of  this  method  use  the  same 
parameters  to  define  a  baseline  problem  as  in  section  3:  M  =  4.0  kg,  m\  =  0.4  kg,  m2  =  1.2  kg, 
L  =  0.254 m,p=  1200  kg/m3, Re(£)  =  108  +105/  N/m2,  andlm(£)=  107  +104/  N/m2, 
where  /is  frequency  in  Hertz.  The  one  additional  parameter  required  is  the  radius  of  the  rotary 
inertial  masses,  which  is  chosen  to  be  0.1016  m.  From  the  above  baseline  values,  the  rotary 
inertial  values  of  the  masses  are  calculated  to  be  J\  =  0.0021  kg»m2  and  /2  =  0.0062  kg»m2, 
which  are  used  to  compute  the  ratios  Aj  and  A2  as  0.247  and  0.741,  respectively.  The  shear 
modulus  values  in  this  analysis  are  Re(G)  =  3.58  x  107  +  3.43  x  104/  N/m2  and  Im(£)  = 

2.55  x  106  +  2.34  x  10 3/  N/m2 ,  where /is  frequency  in  Hertz. 

Figure  10,  corresponding  to  equations  (36)  and  (37),  shows  the  transfer  functions  Sj(£tf)  and 
S2(tt>)  versus  frequency  for  magnitude  (top  plot)  and  phase  angle  (bottom  plot).  The  first 
transfer  function  was  computed  with  an  attached  mass  of  0.4  kg  and  is  depicted  with  x’s;  the 
second  transfer  function  was  computed  with  an  attached  mass  of  1.2  kg  and  is  shown  with  o’s. 

Figure  11  corresponds  to  equation  (40)  and  is  a  plot  of  the  function  r  versus  frequency.  Note 
that  although  this  function  is  a  cosine  with  respect  to  frequency,  the  period  is  increasing  as 
frequency  increases.  Table  3  provides  a  list  of  the  values  of  m  versus  frequency,  which  were 
determined  by  inspection  of  this  figure.  Once  the  values  of  m  are  known,  the  modulus  values  of 
G  can  be  determined  from  equations  (23)-(42). 

Figure  12  is  a  plot  of  real  and  imaginary  shear  modulus  values  versus  frequency.  The  real 
(actual)  values  used  to  make  the  transfer  functions  are  displayed  as  a  solid  line,  and  the  real 
(estimated)  values  are  displayed  as  x’s.  The  imaginary  (actual)  values  used  to  make  the  transfer 
functions  are  also  displayed  as  a  solid  line,  with  the  imaginary  (estimated)  values  displayed  as 
o’s.  The  estimated  values  agree  at  all  frequencies  with  the  actual  values.  This  result  is  expected 
because  there  is  no  noise  in  the  data  and  the  same  parameters  used  to  make  the  transfer  functions 
are  also  used  to  calculate  the  modulus  values  (i.e.,  no  error  is  introduced  when  calculating  the 
modulus  from  the  transfer  functions). 
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Phase 


Frequency  (Hz) 

Figure  10.  Magnitude  and  Phase  Angle  for  Transfer  Functions  Si  and  S2  Versus  Frequency 
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Table  3.  Value  of  m  Versus  Frequency 


Minimum 

Maximum 

Frequency 

Frequency 

(Hz) 

(Hz) 

0 

185 

185 

405 

405 

655 

655 

945 

945 

1280 

1280 

1655 

1655 

2075 

2075 

2540 

2540 

3050 

3050 

3610 

3610 

4215 

4215 

4865 

4865 

5000 

3i 


S] 


Errors,  however,  are  introduced  into  this  method  during  the  measurement  of  the  transfer 
functions  and  system  parameters.  This  effect  can  be  simulated  by  adding  noise  to  the  data  using 
random  numbers  added  to  the  transfer  function,  as  shown  by 

5g(tw)  =  5(£0)  +  e{Re[5(tu)]a,a  +iIm[5(to)<rfe]}  ,  (43) 

where  e  is  the  amount  of  error  added  to  the  transfer  function  and  where  oa  and  Gy  are  random 
numbers  with  zero  mean  and  a  variance  of  one. 

Figure  13,  corresponding  to  equations  (36),  (37),  and  (43),  is  a  plot  of  the  transfer  functions 
Sj(<y)  and  52(<w)  with  noise  versus  frequency.  The  value  of  e  used  to  create  this  figure  is  0.03, 
which  represents  roughly  3%  noise  in  the  data  set.  (If  this  test  is  properly  set  up  and  executed,  a 
3%  noise  level  would  be  considered  greater  than  the  actual  noise  levels  in  the  experimental  data.) 

Figure  14  is  a  plot  of  the  function  r  with  and  without  additive  noise  versus  frequency.  The 
solid  line  is  the  function  with  e  set  equal  to  0.03,  and  the  o’s  connected  with  the  dotted  line  show 
the  function  without  noise  (also  seen  in  figure  11). 

Figure  15  is  a  plot  of  the  estimated  and  actual  values  of  the  shear  modulus  ( G )  versus 
frequency.  The  estimated  values  of  the  real  part  of  G  are  depicted  with  x’s,  and  the  actual  values 
of  the  real  part  of  G  are  shown  as  a  solid  line.  The  estimated  values  of  the  imaginary  part  of  G 
are  depicted  with  o’s,  and  the  actual  values  of  the  imaginary  part  of  G  are  shown  as  a  solid  line. 

The  effect  of  the  measurement  error  on  the  accuracy  of  the  calculation  of  G  can  also  be 
studied  with  Monte  Carlo  simulations.  Eleven  different  values  of  e  were  used  to  build  the 
transfer  functions  and  then  calculate  the  modulus  value.  Estimation  error  at  each  frequency  was 
defined  according  to 


,  .  Gact  )|  Y*est  )| 

a((Om)  = - p - 1 — - Uj  , 

max[GflCf  (Q)m  ),Gest  (o)m )] 


(44) 
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where  cc((Om )  is  the  shear  modulus  estimation  error  at  the  mth  frequency,  the  subscript  act 

corresponds  to  the  actual  value,  and  the  subscript  est  corresponds  to  the  estimated  (calculated) 
value.  Once  this  value  is  known,  it  can  be  summed  across  M  frequencies  by 

1  M 

(45> 

m= 1 

where  is  the  average  error  for  the  y'th  Monte  Carlo  simulation.  Finally,  all  the  simulations  can 
be  summed  from 

e  =  T  X  fy  ’  (4^) 

7=1 

where  e  is  the  estimation  error  using  7  simulations.  Table  4  lists  transfer  function  error  (e)  versus 
shear  modulus  estimation  error  (s).  The  frequency  range  of  the  values  used  to  prepare  this  table 
is  2.5  to  5000  Hz;  the  number  of  frequencies  (M)  used  to  calculate  the  estimation  error  is  100; 
and  the  number  of  simulations  (7)  for  each  value  of  e  is  100. 

Figure  16  is  a  plot  of  shear  modulus  estimation  error  (e)  versus  frequency  for  a  single 
simulation  (7  =  1)  using  a  known  transfer  function  error  of  e  =  0.03.  Note  that  the  estimation 
error  for  the  shear  modulus  is  slightly  larger  than  the  transfer  function  error  in  table  4.  Also  note 
that  the  larger  estimation  errors  tend  to  be  near  the  frequency  values  where  m  is  increasing. 
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Transfer 

Estimation 

Function  Error 

Error 

00 

(e) 

0 

0 

0.005 

0.007 

0.010 

0.013 

0.015 

0.020 

0.020 

0.028 

0.025 

0.035 

0.030 

0.041 

0.035 

0.048 

0.040 

0.052 

0.045 

0.058 

0.050 


0.062 


6.  MEASUREMENT  OF  POISSON’S  RATIO 


The  estimation  of  Poisson’s  ratio  is  achieved  by  combining  the  Young’s  modulus  and  shear 
modulus  previously  measured.  This  equation  is  expressed  as 


v  = 


(47) 


where  v  is  Poisson’s  ratio  and  is  dimensionless.  The  formulation  in  this  method  allows  for 
Poisson’s  ratio  to  be  a  complex  number,  although  typically  the  imaginary  part  of  this  number  is 
very  small  or  zero. 

Figure  17  is  a  plot  of  the  estimated  and  actual  values  of  Poisson’s  ratio  versus  frequency  for  a 
simulation  with  no  noise.  The  estimated  values  of  the  real  part  of  Poisson’s  ratio  are  depicted 
with  x’s,  and  the  actual  values  of  the  real  part  of  Poisson’s  ratio  are  shown  as  a  solid  line.  The 
estimated  values  of  the  imaginary  part  of  Poisson’s  ratio  are  depicted  with  o’s,  and  the  actual 
values  of  the  imaginary  part  of  Poisson’s  ratio  are  shown  as  a  solid  line.  The  estimated  values 
agree  at  all  frequencies  with  the  actual  values,  which  is  expected  because  there  is  no  noise  in  the 
data  and  the  same  parameters  used  to  make  the  transfer  functions  are  used  to  calculate  the  moduli 
values  (i.e.,  no  error  is  introduced  when  calculating  the  moduli  from  the  transfer  functions). 

Figure  18  is  a  plot  of  the  estimated  and  actual  values  of  Poisson’s  ratio  versus  frequency  with 
noise  simulated  using  error  e  =  0.03.  The  estimated  values  of  the  real  part  of  Poisson’s  ratio  are 
depicted  with  x’s,  and  the  actual  values  of  the  real  part  of  Poisson’s  ratio  are  shown  as  a  solid 
line.  The  estimated  values  of  the  imaginary  part  of  Poisson’s  ratio  are  depicted  with  o’s,  and  the 
actual  values  of  the  imaginary  part  of  Poisson’s  ratio  are  shown  as  a  solid  line. 

The  effect  of  measurement  error  on  the  accuracy  of  the  calculation  of  v  can  be  studied  with 
Monte  Carlo  simulations.  Eleven  different  values  of  e  were  used  to  build  the  transfer  functions 
and  then  calculate  the  Poisson’s  ratio.  The  estimation  error  at  each  frequency  was  defined  using 
the  equation 
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(48) 


\Y>act(<0m)\-\Vest(®m)\ 

'X\yJm)  —  r  /  x  T  ? 

max[vact(0)m),vest((0m)\ 

where  a(a)m )  is  the  Poisson’s  ratio  estimation  error  at  the  mth  frequency,  the  subscript  act 
corresponds  to  the  actual  value,  and  the  subscript  est  corresponds  to  the  estimated  (calculated) 
value.  Once  this  value  of  a  is  known,  it  can  be  summed  across  M  frequencies  by 

j  M 

Pj  =  77  X  >  (49) 

1  m=l 


where  fij  is  the  average  error  for  the  jth  Monte  Carlo  simulation.  Finally,  all  the  simulations 
can  be  summed  from 


(50) 


where  e  is  the  estimation  error  using  J  simulations.  Table  5  lists  transfer  function  error  (e)  versus 
Poisson’s  ratio  estimation  error  (s).  The  frequency  range  of  the  values  used  to  prepare  this  table 
is  2.5  to  5000  Hz;  the  number  of  frequencies  (M)  used  to  calculate  the  estimation  error  is  100; 
and  the  number  of  simulations  (J)  for  each  value  of  e  is  100. 


Figure  19  is  a  plot  of  Poisson’s  ratio  estimation  error  (e)  versus  frequency  for  a  single 
simulation  (J  =  1)  using  a  known  transfer  function  error  of  e  =  0.03.  The  estimation  error,  while 
small,  is  not  as  small  as  the  estimation  error  found  in  the  calculation  of  Young’s  and  shear 
moduli.  This  result  is  expected  because  the  errors  from  both  the  estimation  of  Young’s  modulus 
and  shear  modulus  are  contributing  to  the  error  in  the  estimation  of  Poisson’s  ratio. 
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Transfer 

Estimation 
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7.  CONCLUSIONS  AND  RECOMMENDATIONS 


This  report  has  presented  an  innovative  method  for  measuring  the  complex  Young’s 
modulus,  complex  shear  modulus,  and  complex  Poisson’s  ratio  of  a  viscoelastic  cylinder.  The 
new  nonresonant  technique  is  based  on  measured  transfer  functions  that  are  obtained  by  vibrating 
the  cylinder  linearly  and  rotationally  with  two  different-size  masses  on  its  free  end.  Both  masses 
have  their  own  individual  transfer  functions,  which  can  be  measured  on  a  spectrum  analyzer  and 
combined  to  yield  the  unknown  Young’s  modulus  and  shear  modulus  values  at  every  frequency 
where  a  measurement  is  made.  Once  these  moduli  are  determined,  Poisson’s  ratio  can  be 
calculated. 

This  new  method  is  particularly  beneficial  because  it  measures  Young’s  modulus  and  the 
shear  modulus  at  every  frequency  where  a  transfer  function  measurement  is  made  without  relying 
on  system  resonances  or  curve  fitting  to  transfer  functions.  Rather,  the  calculation  from  transfer 
function  measurement  to  the  calculation  of  moduli  is  exact  (i.e.,  no  error  is  introduced). 
Moreover,  Monte  Carlo  numerical  simulations  show  that  this  method  is  relatively  unaffected  by 
noise  introduced  during  the  transfer  function  measurement. 

It  is  suggested  that  future  work  include  laboratory  experiments  to  verify  that  the  closed-form, 
nonresonant  method  can  be  applied  to  actual  physical  problems. 
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